#include <iostream>
#include <fstream>
#include <sstream>
#include <cmath>
#include <random>
#include <exception>
#include "matrix.h"
using namespace std;

static const double min_pos_double = 1e-9;

Shape & Shape::operator=(const Shape & shape){
    rows = shape.rows;
    cols = shape.cols;
    return *this;
}

ostream & operator<<(ostream & os, const Shape & shape){
    os << "(" << shape.rows << ", " << shape.cols << ")" << endl;
    return os;
}

void Matrix::init(int rows, int cols, double ** data){
    rows_ = rows;
    cols_ = cols;
    if(rows == 0 || cols == 0){
        data_ = nullptr;
        return;
    }
    data_ = new double*[rows];
    for(int i = 0; i < rows_; i++){
        data_[i] = new double[cols];
        for(int j = 0; j < cols_; j++){
            data_[i][j] = data[i][j];;
        }
    }
}

bool checkSizeImcompatible(const Matrix & m1, const Matrix & m2) throw(SizeException){
    int rows = m1.getRows();
    int cols = m1.getCols();
    if(rows != m2.getRows() || cols != m2.getCols()){    
        throw SizeException();
        return false;
    }
    return true;
}

bool MatOpSizeImcompatible(const Matrix & m1, const Matrix & m2) throw(SizeException){
    int cols = m1.getCols();
    if(m2.getRows() != cols){    
        throw SizeException();
        return false;
    }
    return true;
}

Matrix::Matrix(int rows, int cols, double default_value){
    rows_ = rows;
    cols_ = cols;
    if(rows == 0 || cols == 0){
        data_ = nullptr;
        return;
    }
    data_ = new double*[rows];
    for(int i = 0; i < rows; i++){
        data_[i] = new double[cols];
    }
    for(int i = 0; i < rows_; i++){
        for(int j = 0; j < cols_; j++){
            data_[i][j] = default_value;
        }
    }
}

Matrix::Matrix(const Matrix & matrix){
    init(matrix.rows_, matrix.cols_, matrix.data_);
}

Matrix::Matrix(Matrix && matrix) noexcept{
    rows_ = matrix.rows_;
    cols_ = matrix.cols_;
    data_ = matrix.data_;
    matrix.data_ = nullptr;
}

Matrix::Matrix(const vector<vector<double>>& value){
    rows_ = value.size();
    cols_ = value[0].size();
    if(rows_ == 0 || cols_ == 0){
        data_ = nullptr;
        return;
    }
    data_ = new double*[rows_]; 
    for(int i = 0; i < rows_; i++){
        data_[i] = new double[cols_];
        for(int j = 0; j < cols_; j++){
            data_[i][j] = value[i][j];
        }
    }
}

Matrix::~Matrix(){
    if(data_ == nullptr) return;
    for(int i = 0; i < rows_; i++){
        delete[] data_[i];
    }
    delete[] data_;
}

double* Matrix::operator[](int idx) const{
    return data_[idx];
    
}

Matrix& Matrix::operator=(const Matrix& matrix){
    if(&matrix == this)
        return *this;
    //free old data
    for(int i = 0; i < rows_; i++){
        delete[] data_[i];
    }
    delete[] data_;
    //malloc new data
    rows_ = matrix.rows_;
    cols_ = matrix.cols_;
    data_ = new double*[rows_];
    for(int i = 0; i < rows_; i++){
        data_[i] = new double[cols_];
        for(int j = 0; j < cols_; j++){
            data_[i][j] = matrix.data_[i][j];
        }
    }
    return *this;
}

Matrix ones_like(const Matrix& matrix){
    return Matrix(matrix.rows_, matrix.cols_, 1);
}

Matrix zeros_like(const Matrix& matrix){
    return Matrix(matrix.rows_, matrix.cols_, 0);
}

ostream& operator<<(ostream& os, const Matrix& matrix){
    for(int i = 0; i < matrix.rows_; i++){
        for(int j = 0; j < matrix.cols_; j++){
            os << matrix.data_[i][j] << TAB_SEP;
        }
        os << LINE_SEP;
    }
    return os;
}

//ostream& operator<<(ostream& os, Matrix&& matrix){
//    for(int i = 0; i < matrix.rows_; i++){
//        for(int j = 0; j < matrix.cols_; j++){
//            os << matrix.data_[i][j] << TAB_SEP;
//        }
//        os << LINE_SEP;
//    }
//    return os;
//}

Matrix & operator+=(Matrix & m1, const Matrix & m2){
    checkSizeImcompatible(m1, m2);
    for(int i = 0; i < m1.rows_; i++){
        for(int j = 0; j < m1.cols_; j++){
            m1.data_[i][j] += m2.data_[i][j];
        }
    }
    return m1;
}

Matrix operator+(const Matrix & m1, const Matrix & m2){
    checkSizeImcompatible(m1, m2);
    Matrix ret;
    double** data = new double*[m1.rows_];
    for(int i = 0; i < m1.rows_; i++){
        data[i] = new double[m1.cols_];
        for(int j = 0; j < m1.cols_; j++){
            data[i][j] = m1.data_[i][j] + m2.data_[i][j];
        }
    }
    ret.init(m1.rows_, m1.cols_, data);
    return ret;
}

Matrix operator+(const Matrix & matrix, const double & num){
    Matrix ret;
    double** data = new double*[matrix.rows_];
    for(int i = 0; i < matrix.rows_; i++){
        data[i] = new double[matrix.cols_];
        for(int j = 0; j < matrix.cols_; j++){
            data[i][j] = matrix.data_[i][j] + num;
        }
    }
    ret.init(matrix.rows_, matrix.cols_, data);
    return ret;
}

Matrix operator+(const double & num, const Matrix & matrix){
    return operator+(matrix, num);
}

Matrix operator-(const Matrix & m1, const Matrix & m2){
    checkSizeImcompatible(m1, m2);
    Matrix ret;
    double** data = new double*[m1.rows_];
    for(int i = 0; i < m1.rows_; i++){
        data[i] = new double[m1.cols_];
        for(int j = 0; j < m1.cols_; j++){
            data[i][j] = m1.data_[i][j] - m2.data_[i][j];
        }
    }
    ret.init(m1.rows_, m1.cols_, data);
    return ret;
}

Matrix operator-(const Matrix & matrix, const double & num){
    Matrix ret;
    double** data = new double*[matrix.rows_];
    for(int i = 0; i < matrix.rows_; i++){
        data[i] = new double[matrix.cols_];
        for(int j = 0; j < matrix.cols_; j++){
            data[i][j] = matrix.data_[i][j] - num;
        }
    }
    ret.init(matrix.rows_, matrix.cols_, data);
    return ret;
}

Matrix operator-(const double & num, const Matrix & matrix){
    Matrix ret;
    double** data = new double*[matrix.rows_];
    for(int i = 0; i < matrix.rows_; i++){
        data[i] = new double[matrix.cols_];
        for(int j = 0; j < matrix.cols_; j++){
            data[i][j] = num - matrix.data_[i][j];
        }
    }
    ret.init(matrix.rows_, matrix.cols_, data);
    return ret;
}
Matrix operator*(const Matrix & m1, const Matrix & m2){
    checkSizeImcompatible(m1, m2);
    Matrix ret;
    double** data = new double*[m1.rows_];
    for(int i = 0; i < m1.rows_; i++){
        data[i] = new double[m1.cols_];
        for(int j = 0; j < m1.cols_; j++){
            data[i][j] = m1.data_[i][j] * m2.data_[i][j];
        }
    }
    ret.init(m1.rows_, m1.cols_, data);
    return ret;
}


Matrix pow(const Matrix& m1, const Matrix& m2){
    checkSizeImcompatible(m1, m2);
    Matrix ret;
    double** data = new double*[m1.rows_];
    for(int i = 0; i < m1.rows_; i++){
        data[i] = new double[m1.cols_];
        for(int j = 0; j < m1.cols_; j++){
            data[i][j] = pow(m1.data_[i][j], m2.data_[i][j]);
        }
    }
    ret.init(m1.rows_, m1.cols_, data);
    return ret;
}

Matrix matmul(const Matrix& m1, const Matrix& m2){
    MatOpSizeImcompatible(m1, m2);
    Matrix ret;
    double** data = new double*[m1.rows_];
    for(int i = 0; i < m1.rows_; i++){
        data[i] = new double[m2.cols_];
        for(int j = 0; j < m2.cols_; j++){
            double element = 0.0;
            for(int m = 0; m < m1.cols_; m++){
                element += m1.data_[i][m] * m2.data_[m][j];    
            }
            data[i][j] = element;
        }
    }
    ret.init(m1.rows_, m2.cols_, data);
    return ret;
}

Matrix operator*(const Matrix & matrix, const double & num){
    Matrix ret = zeros_like(matrix);
    for(int i = 0; i < matrix.getRows(); i++){
        for(int j = 0; j < matrix.getCols(); j++){
            ret[i][j] = matrix.data_[i][j] * num;
        }
    }
    return ret;
}

Matrix operator*(const double & num, const Matrix & matrix){
    return operator*(matrix, num);
}

Matrix Matrix::transpose(){
    Matrix ret(cols_, rows_);
    for(int i = 0; i < cols_; i++){
        for(int j = 0; j < rows_; j++){
            ret[i][j] = data_[j][i];
        }
    }   
    return ret;
}

Matrix exp(const Matrix & matrix){
    Matrix ret(matrix.rows_, matrix.cols_);
    for(int i = 0; i < matrix.rows_; i++){
        for(int j = 0; j < matrix.cols_; j++){
            ret[i][j] = exp(matrix.data_[i][j]);
        }
    }
    return ret;
}

Matrix operator/(const Matrix & m1, const Matrix & m2) throw(DivideZeroException){
    checkSizeImcompatible(m1, m2);
    Matrix ret(m1.getRows(), m1.getCols());
    for(int i = 0; i < m1.getRows(); i++){
        for(int j = 0; j < m1.getCols(); j++){
            ret.data_[i][j] = m1.data_[i][j] / m2.data_[i][j];
        }
    }
    return ret;
}

Matrix operator/(const Matrix & matrix, const double & num) throw(DivideZeroException){
    Matrix ret(matrix.getRows(), matrix.getCols());
    for(int i = 0; i < matrix.getRows(); i++){
        for(int j = 0; j < matrix.getCols(); j++){
            ret.data_[i][j] = matrix.data_[i][j] / num;
        }
    }
    return ret;
}

Matrix operator/(const double & num, const Matrix & matrix) throw(DivideZeroException){
    Matrix ret(matrix.getRows(), matrix.getCols());
    for(int i = 0; i < matrix.getRows(); i++){
        for(int j = 0; j < matrix.getCols(); j++){
            ret.data_[i][j] = num / matrix.data_[i][j];
        }
    }
    return ret;
}

Matrix log(const Matrix & matrix){
    Matrix ret(matrix.rows_, matrix.cols_);
    for(int i = 0; i < matrix.rows_; i++){
        for(int j = 0; j < matrix.cols_; j++){
            ret[i][j] = log(matrix.data_[i][j]);
        }
    }
    return ret;
}

double sum(const Matrix & mat){
    double sum = 0.0;
    for(int i = 0; i < mat.getRows(); i++){
        for(int j = 0; j < mat.getCols(); j++){
            sum += mat[i][j];
        }
    }
    return sum;
}

Matrix readDigitalCSV(string filename){
    ifstream inFile(filename, ios::in);  
    string lineStr;  
    vector<vector<double>> array;  
    while (getline(inFile, lineStr))  
    {  
        // 存成二维表结构  
        std::stringstream ss;  
        ss << lineStr;
        string str;
        vector<double> lineArray;  
        // 按照逗号分隔  
        while (getline(ss, str, ','))  
            lineArray.push_back(atof(str.c_str()));  
        array.push_back(lineArray);  
    }  
    return Matrix(array);
}

Matrix random_matrix(int rows, int cols) {
    vector<vector<double>> ret;
    default_random_engine e;
    uniform_real_distribution<double> u(-1, 1);
    for(int i = 0; i < rows; i++) {
        vector<double> tmp;
        for(int j = 0; j < cols; j++) {
            tmp.push_back(u(e));
        }
        ret.emplace_back(tmp);
    }
    return Matrix(ret);
}
